package odesimulation;

import java.awt.GridLayout;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashMap;
import java.util.HashSet;
import java.util.Iterator;

import javax.swing.JFrame;
import javax.swing.JPanel;

import odesimulation.crnParser.ReactionsetContext;

import org.antlr.v4.runtime.ANTLRInputStream;
import org.antlr.v4.runtime.CommonTokenStream;
import org.apache.commons.math3.ode.FirstOrderDifferentialEquations;
import org.apache.commons.math3.ode.nonstiff.ClassicalRungeKuttaIntegrator;
import org.apache.commons.math3.ode.nonstiff.GraggBulirschStoerIntegrator;

import fitness.ODESineFitness;

import utils.MyStepHandler;
import utils.PlotExpData;

public class OdeImplicitSolver implements FirstOrderDifferentialEquations{

	public static double numberOfSteps = 1000;
	/**
	 * Like OdeSolver, but instead of using an explicit ode system, we only use
	 * an ArrayList of symbols and an ArrayList of Reactions on those.
	 */
	
	ArrayList<String> symbols;
	HashMap<String,Enzyme> enzymes; //map an enzyme to the list of stuff saturating it
	ArrayList<Reaction> reactions;
	
	public OdeImplicitSolver(ArrayList<String> symbols, HashMap<String,Enzyme> enzymes, ArrayList<Reaction> reactions){
		this.symbols = symbols;
		this.enzymes = enzymes;
		this.reactions = reactions;
	}
	
	
	@Override
	public void computeDerivatives(double t, double[] y, double[] ydot) {
		
		Arrays.fill(ydot, 0); //just to be safe, put everything to 0
		HashMap<String, Double> conc = new HashMap<String, Double>();
		for(int i=0; i<y.length; i++){
			conc.put(symbols.get(i), y[i]);
		}
		//Update enzymes
		HashMap<String, HashMap<String, Double>> enz = new HashMap<String, HashMap<String, Double>>();
		for(String e :this.enzymes.keySet()){
			double partialSaturationValue = 1.0;//Still has to be devided in the final value by the relevant Km
			Enzyme actualEnzyme = this.enzymes.get(e);
			HashMap<String, Double> actualKis = new HashMap<String, Double>();
			for(String sat : actualEnzyme.Kms.keySet()){
				partialSaturationValue += conc.get(sat)/actualEnzyme.Kms.get(sat);
			}
			for(String sat : actualEnzyme.Kms.keySet()){
				actualKis.put(sat, actualEnzyme.kis.get(sat)/partialSaturationValue);
			}
			enz.put(e, actualKis);
		}
		
		for(Reaction r : reactions){
			HashSet<String> relevantSpecies = new HashSet(r.products); //so we remove double
			relevantSpecies.addAll(r.reactants);
			Iterator<String> it = relevantSpecies.iterator();
			while(it.hasNext()){
				String s = it.next();
				if(r.getClass() != EnzymaticReaction.class){
				ydot[symbols.indexOf(s)] += r.computeFlow(s, conc);
				} else {
					String enzyme = ((EnzymaticReaction) r).enzyme;
					Enzyme actualEnzyme = this.enzymes.get(enzyme);
					String origin = r.reactants.get(0);
					double RelevantKm = actualEnzyme.Kms.get(origin);
					double RelevantKi = enz.get(enzyme).get(origin);
					
					ydot[symbols.indexOf(s)] += ((EnzymaticReaction) r).computeFlow(s, RelevantKi/RelevantKm, conc);
				}
			}
		}
	}

	@Override
	public int getDimension() {
		return symbols.size();
	}

	public static ArrayList<Reaction> generateReactions(ReactionsetContext context){
		ArrayList<Reaction> reaclist = new ArrayList<Reaction>();
		
		return reaclist;
	}
	
	/**
	 * Solve a specific set of differential equations
	 * @param initialConditions Initial concentrations for all substrats
	 * @param symbols Names of substarts
	 * @param enzymes Enzymes present in the system. Link between their name and their physical representation
	 * @param reactions All the reactions present in the system.
	 * @return
	 */
	public static double[][] solveEqus(double[] initialConditions, ArrayList<String> symbols,  HashMap<String,Enzyme> enzymes, ArrayList<Reaction> reactions){
		
		OdeImplicitSolver odes = new OdeImplicitSolver(symbols, enzymes, reactions);
		GraggBulirschStoerIntegrator myIntegrator = new GraggBulirschStoerIntegrator(
				1e-13, 1, 1e-6, 1e-6);
		MyStepHandler handler = new MyStepHandler(odes.getDimension());
		myIntegrator.addStepHandler(handler);
		myIntegrator.addEventHandler(new utils.CalcErrEventHandler(), 2000, 0.000000001, 100000);
		myIntegrator.addEventHandler(new utils.MyODEEventHandler(), 10, 1, 100);
		myIntegrator.integrate(odes, 0, initialConditions, OdeImplicitSolver.numberOfSteps ,
				initialConditions);
		return handler.getTimeSerie();
	}
	
	public static double fitness1(String s, int[][] target){
		double value = 0.0;
		int sizex = target[0].length;
		int sizey = target.length;
		
		//Lexing the equations
		crnLexer myLex = new crnLexer(new ANTLRInputStream(s));
		//Parsing the equations
		crnParser myPars = new crnParser(new CommonTokenStream(myLex));
		//Parsing result
		//Contains "symbols": X,Y,Z... list of all substrats.
		//		   "enzymes": E1... list of all enzymes.
		//         "reactions": list of all reactions.
		ReactionsetContext context = myPars.reactionset();
		
		//initial concentrations for substrats: 10 nM all around!
		double[] init = new double[context.symbols.size()];
		Arrays.fill(init, 5);
		
		return value;
	}
	
	public static void main(String[] args){
		//Lexing the equations
		double k1 = 0.1;
		double maxk1 = 4.0;
		int sizex = 5;
		int sizey = 5;
		double mink4 = 0.0;
		double maxk4 = 8.0;
		double stepx = (maxk4-mink4)/((double)sizex);
		double stepy = (maxk1-k1)/((double)sizey);
		
		String exp = "";
		JFrame frame = new JFrame();
		JPanel content = new JPanel(new GridLayout(0,sizex));
		double time = System.currentTimeMillis();
		while(k1 < maxk1){
		double k4 = mink4;
			while(k4 < maxk4){
				
		crnLexer myLex = new crnLexer(new ANTLRInputStream(
				"# individual \"99\"\n"+
				"enzymes: E1 E2 E3"+
				"products: X Y Z"+
				"K_m[3][3] : { { 35.35651646265842 1.0 18.074538782097136 }{ "+(1.0/(k1*8.7e-2))+" 0.19639942990693068 1118.6216070905916 }{ 49.06140644009669 7.5529597327930205 4.656470904572082 } }"+
				"k_i[3][3] : { { 5392.936282100518 0.39 3.3851197232870134 }{ "+(1.0/8.7e-2)+" 6.127154985599243 7.02027575785981 }{ 275.49850207203826 16.44430374780971 2.3695397213682563 } }"+
				"X + E2 ==> E2 + X+X "+
				"X + Y <==> Y+Y [ k+ = 1.0 ; k- = 0.0 ]"+
				"X <==>   [ k+ = "+0.39*k4+" ; k- = 0.000 ]"+
				"Y + E1 ==>  E1"+
				" <==> X [ k+ = 0.0000001 ; k- = 0.0 ]"));
		//Parsing the equations
		crnParser myPars = new crnParser(new CommonTokenStream(myLex));
		//Parsing result
		//Contains "symbols": X,Y,Z... list of all substrats.
		//		   "enzymes": E1... list of all enzymes.
		//         "reactions": list of all reactions.
		ReactionsetContext context = myPars.reactionset();
		
		//initial concentrations for substrats: 10 nM all around!
		double[] init = new double[context.symbols.size()];
		Arrays.fill(init, 5);
		
		double[][] result = solveEqus(init,context.symbols,context.enzymes,context.reactions);
		//System.out.println(k1+" "+k4);
		fitness.OscillatorBasedResult fit = ODESineFitness.evaluateFitness(result, context.symbols);
		exp += (fit.finalResult()>1?" "+fit.average+"\t":" -\t");
		String[] seqNames = new String[context.symbols.size()];
		seqNames =context.symbols.toArray(seqNames);
		
		PlotExpData demo = new PlotExpData("",
				result, seqNames, "", "");
		
		content.add(demo);
		
		k4 += stepx;
			}
			//System.out.println(exp);
			exp += "\n";
			k1+= stepy;
		}
		double time2 = System.currentTimeMillis();
		System.out.println(exp);
		System.out.println("\n"+(time2-time));
		frame.add(content);
		frame.pack();
		frame.setVisible(true);
	}
}
